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ABSTRACT 

We present a simple and powerful method for extracting transit signals associated with a known transiting 
planet from noisy light curves. Assuming the orbital period of the planet is known and the signal is periodic, 
we illustrate that systematic noise can be removed in Eourier space at all frequencies, by only using data within 
a fixed time frame with a width equal to an integer number of orbital periods. This results in a reconstruction 
of the full transit signal which on average is unbiased, despite that no prior knowledge of either the noise or 
the transit signal itself is used in the analysis. The method has therefore clear advantages over standard phase 
folding, which normally requires external input such as nearby stars or noise models for removing systematic 
components. In addition, we can extract the full orbital transit signal (360 degrees) simultaneously, and Kepler 
like data can be analyzed in just a few seconds. We illustrate the performance of our method by applying it to 
a dataset composed of light curves from Kepler with a fake injected signal emulating a planet with rings. Eor 
extracting periodic transit signals, our presented method is in general the optimal and least biased estimator 
and could therefore lead the way toward the first detections of, e.g., planet rings and exo-trojan asteroids. 


1. INTRODUCTION 

Our understanding of where exoplanets form, their or¬ 
bital distribution, and even properties such as mass and ra¬ 
dius has greatly improved over recent years (Lissauer et al. 
2014; Laughlin & Lissauer 2015). Several methods have 
contributed to this progress, from mass sensitive measure¬ 
ments including pulsar timings (Wolszczan & Erail 1992), ra¬ 
dial velocimetry (Mayor &. Queloz 1995; Mayor et al. 2014) 
and microlensing (Bond et al. 2004), to transit photome¬ 
try (Henry et al. 2000; Charbonneau et al. 2000), which in 
turn depends on the planet radius. Recent developments 
have primarily been driven by photometric surveys operat¬ 
ing from both the ground (OGLE (Udalski et al. 2002), HAT 
(Bakos et al. 2002) and WASP (Pollacco et al. 2006)) and 
space (CoRoT (Auvergne et al. 2009), Kepler (Borucki et al. 

2010) ). The emerging picture is that exoplanets are not 
rare, but exist with a wide variety of masses, orbital sizes, 
and eccentricities, and often in multiplanet configurations 
which can be very different from our own solar system 
(Laughlin & Lissauer 2015). Despite the success in detecting 
planets, there are still major gaps in our understanding of how 
they form and their subsequent evolution. The problem covers 
all aspects of astrophysics from the formation and evolution 
of the protoplanetary disk (see e.g.. Turner et al. (2014)), to 
late time gravitational dynamics including both strong field 
planet-planet scattering and secular evolution (Beauge et al. 
2012). Due to this complexity, observations play a crucial 
role in driving the field. 

In this paper we focus on photometric measurements and 
how to improve the scientific output from transit light curves. 
One of the present challenges is to reduce the noise to 
a level where ’higher-order photometric effects’ (Perryman 

2011) , with relative flux variations of order 10“^, can be re¬ 
solved. This level of precision could lead to the first de¬ 
tection of exomoons (Simon et al. 2012; Kipping et al. 2012; 
Heller 2014; Hippke 2015), planet rings (Barnes & Eortney 
2004), and exo-trojan asteroids (Laughlin & Chambers 2002), 
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as well as teach us about the rotation of the transiting plan¬ 
ets (Zhu et al. 2014) and even their atmosphere (Hui & Seager 
2002; Sidis & Sari 2010), to name a few examples. However, 
such variations are very hard to resolve with present observa¬ 
tions and noise removal techniques. 

Motivated by these challenges, we here present a fast and 
effective method for extracting small flux variations associ¬ 
ated with a known planet transit, from light curves with ran¬ 
dom and systematic noise. The noise represents here all com¬ 
ponents which are not associated with the transit signal it¬ 
self such as instrumental noise and real stellar flux varia¬ 
tions. By assuming the transit signal is periodic, we illustrate 
how one can remove systematic noise in a model independent 
way, which makes it possible to resolve finer flux variations 
than standard phase folding techniques^. This could directly 
improve the ongoing search for planetary rings and trojan- 
asteroids, as well as lead to better measurements of any phe¬ 
nomena associated with the planet transit itself. 

The new strategy we propose is based on removing system¬ 
atic noise in Eourier space. We show this becomes especially 
effective if one only uses data within a fixed time frame start¬ 
ing from the middle of the first planet transit and ending in 
the last. In this way any signal associated with the transiting 
planet will be constrained to discrete and separate frequen¬ 
cies in Eourier space, which makes the noise separation pos¬ 
sible. Our method does not rely on any prior knowledge of 
either the noise or signal, which also makes it fast enough 
that several quarters of Kepler data can be analyzed in sec¬ 
onds. If the noise is not directly correlated with the periodic 
signal, our method will, on average, return an unbiased recon¬ 
struction of the full orbital transit signal. Eourier transforms 
have previously been used for, e.g., finding short period plan¬ 
ets (Sanchis-Ojeda et al. 2014), as well as removing noise in 
exoplanet emission spectra (Swain etal. 2010; Zellemet al. 
2014), but has to our knowledge not been applied for extract¬ 
ing the transit shape itself in the way we propose. 

This paper is organized as follows. In Section 2 we describe 

^ Phase folding techniques are in this paper loosely refening to methods 
which use out-of-transit fits and stacking in real space to reduce systematic 
and random noise, respectively. 
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how the observed flux, possible signal, and noise are related. 
The notation and assumptions from this section will be used in 
the remaining parts. In Section 3 we describe our new Fourier 
method step by step and how it can be applied to reconstruct 
any periodic transit shape. The performance of the method is 
shown in Section 4, where we apply it to reconstruct a fake 
signal injected into real Kepler light curves. Conclusions are 
given in Section 5. 

2. TRANSIT SIGNALS AND LIGHT CURVE OBSERVATIONS 

We consider a planetary system consisting of a star with a 
single planet. If observed from the right position, the planet 
will move in front of the star during its orbit, which gives rise 
to periodic dips in the observed light curve. We denote this 
the main transit. Objects associated with the planet such as 
rings, moons, and trojan-asteroids can likewise block parts of 
the stellar light, and can therefore be indirectly seen as small 
flux variations on top of the main transit signal. A few ex¬ 
amples are shown in Figure 1. In this section we describe 
the assumptions we make about how the observed flux, transit 
signal, and noise are related to each other. 

2.1. Flux from a Star with an Orbiting Planet 

For this paper we consider observations of a star with one 
transiting planet. The observed flux Fobs can in the majority 
of cases be modeled as 

Fobs=F* + AF + F*/‘S, (1) 

where F,. is the mean flux from the star, AF is considered 
noise, and is the fractional flux modulation coursed by or¬ 
biting objects including the planet. The modulation term 
is typically negative, but reflections from the planet surface 
can give rise to positive modulations. The noise term AF 
is assumed to be small compared to Ff, and represents all 
systematic and random flux variations which are not associ¬ 
ated with the true transit signal. By defining the observable 
quantity I{t) = Fobs/F, - 1 and the relative noise contribution 
/(f)^ = AF/Ff, we can rewrite equation (1) as 

I{t) = I{tf+I{tf. (2) 

In the remaining part of the paper we generally refer to I{t) 
as the data, I(t)^ as the signal and I(t)^ as the noise. In the 
following sections we show that if the signal I{t)^ is periodic, 
then one can make a precise reconstruction of its transit shape 
by working in Fourier space. 

3. RECONSTRUCTION OF PERIODIC TRANSIT SIGNALS 




Fig. 1.— Illustration of a planet with rings and a population of exo-trojan 
asteroids. The top plot shows a schematic drawing of the system, where 
the bottom plot shows the con'esponding transit signal. The illustration is 
not to scale, but intends to illustrate how a typical signal could look like 
and especially how planet rings and asteroids contribute to the signal. A 
real signal has encoded much more information such as how fast the planet 
rotates and if it has an atmosphere, to name a few examples. However, it is in 
general very hard to resolve the resultant flux valuations which would show 
up as small perturbations to the main transit (solid line) at the percent level. 
In this paper we present a way of improving this by working in Fourier space. 


represented as a Fourier series 



( 3 ) 


where t is time, N is the number of data points, and r,, Tf 
are denoting the starting time and ending time of the dataset 
used in the analysis, respectively. We refer to this fixed time 
frame as the ’data window’ and denote its width, |ry -r,j, by 
T„. The importance of the width and location of this window 
is discussed in section 3.2. From our signal and noise model 
shown in equation (2) and by the linearity of the Fourier trans¬ 
formation, the Fourier coefficients can be divided up into sig¬ 
nal and noise as 

an = af + a^, bn = b^ + b^. (4) 

The idea is now to separate the signal and the noise in the 
Fourier domain, as shown in equation (4). In the section be¬ 
low we show that a specific choice for r, and tj will make this 
approach possible. 


We now present our new method for extracting periodic 
transit signals. The method relies on the basic assumption 
that the transit signal must be strictly periodic, with the same 
period as the orbiting planet. Besides that, there are no restric¬ 
tions on the shape of the signal, i.e., the method can be used to 
reconstruct transits with any shape. The assumed periodicity 
makes it possible to work in Fourier space, which enables us 
to remove systematic noise components. This is in contrast to 
standard phase folding techniques where the real space stack¬ 
ing is only optimal for reducing random noise. In the sections 
below we describe our proposed Fourier method step by step. 

3.1. Representing Data as a Fourier Series 

The general strategy is to extract the periodic transit sig¬ 
nal in Fourier space. We therefore first consider the data /(f) 


3.2. Choosing the Optimal Data Window 

The Fourier coefficients obtained from using the full data 
set /(f) will in general have broad features, and the underly¬ 
ing signal will overlap with the noise spectrum, even if the 
signal has strictly periodic components. This is mainly due 
to ’window effects’ where the position and width of the data 
window interfere with the Fourier modes of the data. Separat¬ 
ing the coefficients as illustrated in equation (4), is therefore 
not possible in general. 

However, what is usually not appreciated is that the data 
window actually can be chosen freely to highlight features of 
interest in the data. In our case, we can use our knowledge of 
where the main transits occur in the light curve to construct a 
window which will not affect the representation of the signal 
in Fourier space. Following this idea, one can show that the 
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minimum overlap between noise and signal generally can be 
achieved by choosing a data window such that 

r, = time of first transit 
Tf = time of last transit, 

where the transit time here refers to the time at the center of 
the main planet transit. There are two reasons for why this 
window can significantly improve the reconstruction of the 
signal: first, the frequencies related to the signal will now be 
confined to single frequency bins with no broadening because 
T„ is now a multiple of the period of the planet. Second, since 
the time window is now fixed to start in the middle of a transit, 
the obtained signal coefficients, af and , not only describe 
the odd and even contributions to the data, but also the sym¬ 
metry of the signal itself. For example, if we search for a 
signal which is symmetric around the main transit (as the one 
shown in Figure 1), then we know right away that all the hf 
must equal zero. This gives us a perfect separation of noise 
and signal in the coefficient space, which greatly improves 
the extraction of the signal. The top plot in Figure 2 illustrates 
the proposed window applied to data from Kepler, where the 
top center plot shows how this results in Fourier coefficients 
with distinct and isolated signal peaks. 

3.3. Reconstruction of a Symmetric Transit Signal 

To make the presentation of the method as clear as possi¬ 
ble and to shorten the algebra, we now assume that the tran¬ 
sit shape of the signal is symmetric with respect to the main 
transit. This is a special case, but in fact represents a large 
class of real astrophysical transit phenomena including the 
main planet transit, secondary eclipse, planetary rings, and 
trojan-asteroids. Figure 1 shows an example. In section 3.4 
we describe how to extend the following procedure to recon¬ 
struct a signal with arbitrary transit shape. 

3.3.1. Identifying Noise and Signal in Fourier Space 

By assuming the transit shape is symmetric, the signal peri¬ 
odic, and using the window proposed in equation (5), we then 
know that only cosine terms will be needed for the reconstruc¬ 
tion. Inserting this case in equation (3) and using the linearity 
relation between signal and noise shown in equation (4), we 
can now express the data described by the cosine terms as 



(6) 


where Nq denotes the number of full planet orbits within the 
data window. From this equation it is now clear that the noise 
can contribute at all n = {1,2,3...}, but the signal only con¬ 
tributes at n = mNo where m = {1,2,3...}. In other words, the 
signal is now constrained to discrete and separate frequency 
bins in Fourier space. The fraction these signal bins cover of 
the total available frequency space is ^ 1 /{2Nq), which is sim¬ 
ply the ratio of the number of signal bins, N/(INq), to the total 
number of bins, N. This means that if A^o = 10 then the sig¬ 
nal only appears in 5% of the total available space, i.e., 95% 
of the modes in Fourier space can be labeled as pure noise. 
An example of this signal and noise separation is shown in 
the middle panels of Figure 2. In the next section we illustrate 
how the noise contribution can be estimated in the signal bins, 
which leads to the final extraction of the signal. 


Kepler: 010264660-2009259160929 (REAL Transit Signal) 
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Fig. 2.— Extraction of a planet transit signal from real Kepler data using 
our proposed Fourier method. The steps in the method are shown from top 
to bottom. Top: Rescaled PDCSAP flux from the Kepler observed system 
010264660-2009259160929. The gray line indicates the full dataset from 
that quarter, where the black points are the measurements we select for our 
analysis. As seen, these points are all within the vertical dashed lines which 
are chosen to go through the exact middle of the first and the last planet tran¬ 
sit, respectively. This specific data window greatly improves the separation 
of signal and noise in Fourier space as described in Section 3.2. Top center: 
Fourier coefficients an obtained from representing the selected dataset in a 
Fourier series, as shown in equation (3). The frequency bins with clear nar¬ 
row peaks pointing toward negative values are the signal bins and contain a 
mix of noise and signal. All other bins contain pure noise. Bottom center: 
The Fourier coefficients an separated into a noise component a{^ (red) and a 
signal component af (blue). The noise contribution in the signal bins is here 
estimated by using linear inteipolation as described in Section 3.3.2. This 
way of removing systematic noise makes our method on average an unbiased 
estimator of the true transit signal, in contrast to nomial phase folding. Bot¬ 
tom: The full reconstructed signal shown up to the Nyquist frequency. The 
blue symbols denote a^ (even part) coefficients, where the green symbols 
denote bf (odd part) coefficients. 
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3.3.2. Reconstruct Signal 

From using the fixed data window proposed in section 
3.2, the signal has now been constrained to only appear in 
^ 1 /(2Nq) of the available Fourier space. However, noise will 
in general still be present in the bins where the signal is now 
localized. The question is if we can remove this noise con¬ 
tribution to improve the reconstruction of the true signal. We 
start by using equation (4) to express the signal coefficients as 

S J\f cn\ 

^mNo ~ ^>"Ro ~ ^mNo ’ s') 

where a^No is obtained directly from data. If one now can es¬ 
timate then the underlying signal can be extracted. 
The random part of the noise can of course not be removed 
in individual bins, but we here show that it is possible to cor¬ 
rect for systematic components since these on average appear 
as smooth features in Fourier space. The idea is to use the 
spectral shape of the noise outside the signal bins {n ^ niNo) 
to predict the value inside the bins (n = niNo). In other words, 
we suggest that the systematic contribution to the noise spec¬ 
trum can on average be estimated by interpolating across the 
signal bins in Fourier space, using the regions outside the bins 
as a baseline. The regions over which we have to interpolate 
are always just single frequency bins with unit width, no mat¬ 
ter how broad the signal is in real space. The optimal inter¬ 
polation scheme depends on the properties of the noise, how¬ 
ever, we find that a simple linear interpolation is by far the 
most robust when using real data. Using a linear interpolation 
scheme, the estimate for the true signal coefficients is now 
given by 

‘^mNo ~ ‘^niNo ~ («mVo-l + OmNo+l)/'2', (8) 

where the last term is just the average of the neighbor bins 
around The coefficients can directly be found from 
data by 


N 

ax = y^J{tj)cos 

J=i 



0<tj< 


( 9 ) 


where the summation must be performed over the fixed time 
window as indicated. If the data is nonuniformly sampled 
in time, one must use a slightly more general approach for 
estimating the coefficients. This is further discussed in the 
Appendix. 

It is often useful to also consider the real space representa¬ 
tion of the reconstructed signal. This representation is simply 
given by 



where Tq is the orbital time of the planet, and the signal coeffi¬ 
cients are given by equation (8). The real space represen¬ 
tation has the advantage that the different signal components 
can separate in time along the planetary orbit (see e.g.. Figure 
3). However, the drawback is that single outliers among the 
Fourier coefficients generate global periodic variations, which 
makes it very hard to distinguish noise in Fourier space from 
a possible interesting signal identified in real space. This is 
different from the Fourier space representation where all the 
signal components share the same frequency bins. A model 
comparison must therefore be performed by fitting in Fourier 
space instead of real space. A blind search for unexpected 


Kepler: 010264660-2009350155506 (FAKE Transit Signal) 
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Fig . 3.— Illustration of how well our proposed Fourier method reconstructs 
a fake signal injected into real light curves from the Kepler telescope. For this 
example we use 12 quarters of data from the dataset 010264660. In addition 
to the main transit from the planet, we have also added a slightly broader 
component with a transit depth of 10% relative to the main transit to emulate 
a ring system (Barnes & Foitney 2004; Tusnski &, Valio 2011). The top plot 
shows one of the 12 quarters of data including the fake injected signal. We use 
the data inside the matched window illustrated by the vertical dashed lines, 
as described in Section 3.2. The top center panel shows the corresponding 
Fourier coefficients where a,, is shown in blue and bn in green. The median 
values from all 12 quarters are shown in black and light green, respectively. 
In the right plot is shown the Fourier representation of the injected fake signal 
with a red solid line, where the orange dashed shows the contribution from 
the planet and the orange dashed-dotted the contribution from the rings. The 
bottom center plot shows the real space representation of the reconstructed 
signal using the median estimate. The bottom plot shows a zoom in on the 
main transit. As seen, our method (black symbols) reconstructs the true shape 
(solid red) to a precision where the difference between a signal without rings 
(dashed line) and with rings (solid line) clearly can be seen. The same is true 
for the Fourier representation shown in the top center plot. Any transit shape 
can be reconstructed by this Fourier method, as long as the signal is strictly 
periodic. 
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signals is on the other hand probably best to perform in real 
space, where the different parts of the signal often separate 
and therefore becomes easier to visually identify. 

Our estimate of the signal shown in equation (8) is of course 
not always perfect, but it is expected to be good on average. 
This means that if one have several datasets from the same 
system, the stacked set of reconstructed signals in Fourier 
space will on average be an unbiased estimate of the true sig¬ 
nal. Standard phase folding is not unbiased in the same way, 
since the stacking in real space is only optimal for reducing 
random uncorrelated noise. A noise contribution with, e.g., a 
sinusoidal shape has highly correlated values in real space and 
is therefore not averaging to zero if divided up in orbital seg¬ 
ments and stacked, unless the divisions are matched in phase 
and size to the noise shape. However, in Fourier space such a 
contribution is basically only occupying a single bin and can 
therefore in general be fully removed without touching the 
signal bins. Masking out the signal regions in real space and 
using the data outside to estimate systematic noise is some¬ 
times possible in phase folding, but highly limited to cases 
where the location of the signal is already known and con¬ 
strained to a small region out of the full orbit. These restric¬ 
tions do not apply to our method, which makes it possible to 
remove noise at all frequencies and reconstruct the full orbital 
signal. Our new method is illustrated step by step in Figure 2 
on Kepler data. 

3.4. Reconstruction of an Arbitrary Transit Signal 

Any function can be described by the sum of an even and 
an odd function. In this section we have shown how the even 
part of the transit shape can be reconstructed. The reconstruc¬ 
tion of the odd part is identical to this procedure, except that 
sine terms are then used instead of cosine terms. Any transit 
shape can therefore be extracted from a dataset with both ran¬ 
dom and systematic noise using our proposed method, with 
the only requirement that the signal must be strictly periodic. 

4. NUMERICAL EXAMPLE: A PLANET WITH RINGS 

We here give an example of how well a transit signal can 
be reconstructed using our Fourier method. For this illustra¬ 
tion we inject a fake signal into 12 quarters of data from the 
Kepler object 010264660. This object already has a planet, 
as seen in Figure 2, which we therefore remove before inject¬ 
ing our fake signal. Besides the clean transit signal from the 
planet itself, we also include a wider and weaker transit com¬ 
ponent centered around the planet to emulate a ring system. 
This gives rise to wide wings at ingress and egress, which is 
one of the main indications for rings (Barnes & Fortney 2004; 
Tusnski & Valio 2011). For illustrative purposes, we chose 
the injected signal to have a symmetric shape with respect to 
the main transit. The interesting question is now if our method 
is precise enough to tell if a ring is present or not. 

Results are summarized in Figure 3. There are several ways 


of how to combine the 12 datasets, here we use the approach 
of applying our method to each of the datasets individually. 
As seen in the panel showing the Fourier coefficients, each 
realization has some scatter, but the average of the 12 re¬ 
constructed sets is consistent with the true signal shape. As 
a result, the contribution from rings are easily seen both in 
Fourier space (top center panel) and in real space (bottom 
panel). How to compare the extracted signal to models in¬ 
cluding proper error propagation, will be discussed in an up¬ 
coming paper. 

5. CONCLUSIONS 

We present a simple, yet powerful, approach for extracting 
the shape of periodic transit signals associated with a known 
planet from an observed light curve. The method is based on 
the idea of separating signal and noise in Fourier space, which 
especially enables us to remove systematic noise components. 
We show that this becomes possible if only data within a fixed 
time window with a width equal to an integer number of or¬ 
bital periods, is used in the analysis. In this case, the sig¬ 
nal separates out and is described by a significantly smaller 
number of frequencies than the full dataset. In addition, the 
signal frequencies are easily identified in Fourier space, and 

1 -(1 /Nf) of the total space can be labeled as pure noise as 
a result. For extracting the true transit signal, we suggest that 
the noise contribution across the frequency bins of the signal 
can be estimated by using the regions outside the bins as a 
baseline. We illustrate this by successfully extracting a fake 
transit signal injected into real Kepler light curves. 

Our proposed method for removing systematic noise does 
not rely on either input noise models or nearby field stars for 
calibration. This makes our method not only model indepen¬ 
dent, but also fast enough that several quarters of Kepler data 
can be analyzed in seconds. Furthermore, we can extract the 
full orbital transit signal without applying any masking or out- 
of-transit fitting. Our method has of course some limitations 
which depend on the data quality and the noise properties, 
and one could in general also benefit from combining it with 
other techniques such as the complementary trend filtering al¬ 
gorithm (TFA) (Kovacs et al. 2005). This will be explored in 
an upcoming paper where we apply the method to real data, 
to search for rings, trojans and other higher-order photometric 
effects associated with main orbiting planet. 
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G. Bakos, S. Pedersen, C. Huang, M. Val-Borro, T. Brandt, J. 
Hartman for useful conversations, and especially C. Holcomb 
for going through an earlier version of the manuscript. Sup¬ 
port for this work was provided by NASA through Einstein 
Postdoctoral Fellowship grant number PF4-150127 awarded 
by the Chandra X-ray Center, which is operated by the Smith¬ 
sonian Astrophysical Observatory for NASA under contract 
NAS8-03060. 


APPENDIX 

Observationally obtained light curves are often sampled nonuniformly in time and include gaps with no data. The Fourier 
coefficients can therefore not be obtained using standard FFT methods, nor by the simple relation shown in equation (9). Instead, 
one has to perform a joint least square fit of sine and cosine functions to the data. Due to the linearity of the Fourier expansion, 
this fit is simple and can with great advantage be formulated using matrices. Following this approach we can now express the 
data as a matrix multiplication 

I; K, C„;an-P (1) 

where C„; = cos(27r«f;/r„), S„, = sin(27rnf,/T'„), Tj-,- is their joined matrix, and = [a„,b„] is the vector of coefficients we want to 


6 


fit for. By minimizing the difference squared between our expansion and the data, one can now express the solution to the vector 
Ra- as 

Ra ft! (T'^T)“‘tTI, (2) 

where the superscript T here denotes the transpose. If the data is uniformly sampled then the matrix T^T equals the identity 
matrix, and the solution reduces to the simple Fourier relation shown in equation (9). Only in this case can a FFT solver be used. 
Data from the Kepler telescope is very close to being uniformly sampled, however there are many gaps with missing data. One 
must therefore use this general fitting approach to keep the frequencies under exact control. All examples in this paper are based 
on this general fitting. 
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